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Abstract. We present a detailed investigation of the cosmological dynamics 
based on exp (-R/A) gravity. We apply the dynamical system approach to both 
the vacuum and matter cases and obtain exact solutions and their stability in the 
finite and asymptotic regimes. The results show that cosmic histories exist which 
admit a double de-Sitter phase which could be useful for describing the early and 
the late-time accelerating universe. 



1. Introduction 

In recent years several observational surveys appear to show that the universe is 
currently undergoing an accelerated expansion phase (for a review see pQ). On the 
other hand, astrophysical observations on galactic scales gave us a clear indication 
that the amount of the luminous matter is not enough to account for the rotation 
curves of galaxies. These remarkable observations have radically changed our ideas 
about the evolution of the universe. In fact, the only way in which this behaviour can 
be explained within the standard Friedmann general relativistic cosmology framework 
is to invoke two dark matter components, one with a conventional dust equation of 
state (w = 0) needed to fit astrophysical data on galactic scales (Dark Matter) and 
the other with a more exotic equation of state (iu < — |) needed to explain the current 
accelerated expansion phase of the universe (Dark Energy) . There is little doubt that 
unraveling the nature of Dark matter and Dark Energy is currently one of the most 
important problems in theoretical physics. 

One approach to the problem of Dark Energy that has received considerable 
attention in the last few years is the modification of General Relativity on cosmological 
scales. In this approach one supposes that Dark energy is a manifestation of a non- 
Einsteinian behaviour of the gravitational interaction rather than a new form of energy 
density. The introduction of corrections to the Hilbert Einstein action and their effects 
have been studied for long time and are believed to be unavoidable when the quantum 
nature of the universe is introduced in General Relativity (GR). 

In the last few years many different extended versions of the Einstein theory 
of gravity have been proposed. One of the most studied approaches is higher order 
theories of gravity [13 — 42], in which the gravitational action is nonlinear in the 
Ricci curvature and/or it's derivatives El [8]. These theories have a number of 
interesting features on cosmological and astrophysical scales. In fact they are known 
to admit natural inflation phases [4] and to explain the flattening of the galactic 
rotation curves [9] . Another very interesting feature of these models is that the higher 
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order corrections to Hilbert-Einstein action can be viewed as an effective fluid which 
can mimic the properties of Dark Energy [5] . 

One of the the main problems that occurs in the study of higher order theories 
of gravity is that finding exact cosmological solutions is extremely difficult due to 
the high degree of non-linearity exhibited by these theories. This problem can be 
partially addressed using a suitable choice of generalized coordinates in which the field 
equations can be written as a system of first-order autonomous differential equations 
together with a constraint equation [10] . In this way we can exploit the methods of 
dynamical systems theory [2] in order to both understand the qualitative behavior 
of the cosmological dynamics and obtain special exact solutions of the cosmological 
equations. The general approach allowing one to analyze higher order gravity with 
dynamical systems techniques has been presented elsewhere [3]. 

In this paper we will apply this approach to an important class of theories with 
Lagrangian density 

e-%+L M , (1) 

where A is the cosmological constant. This Lagrangian has some interesting features, 
first of all its treatment is equivalent to one involving a combination of powers of the 
Ricci scalar and, secondly it reduces to 

ex P (-£) = l-|+0[i? 2 ], (2) 

In the small curvature limit, which is equivalent to Hilbert-Einstein action. 

This paper has been arranged as follows. In section 2, we present the basic 
equations of the model. In section 3 and 4 we find exact solutions and their stability in 
the vacuum and matter cases respectively. Finally in section 5 we present a discussion 
of the results and present our conclusions. 

In what follows we will use natural units (H = c = kg = 8ttG = 1) and the 
signature (+,—,—,—). 



A 



d 4 xC 



2. Basic Equations 



The general action for a fourth order theory of gravity is: 

A= ( d 4 x^g[f{R) + L M ], (3) 



where f(R) is a function of Ricci scalar R and Lm is the standard matter Lagrangian 
density. Varying this action with respect to the metric gives 

p _ rpToT _ rpM , rpR _ 1 rpM 

+ -Jm (\g^ lf(R) - Rf\R)] + f\Ry aP tecvgpv - 9^9^) ) , (4) 

where the prime denotes the derivative with respect to R, T^f v is the effective stress- 
energy tensor for standard matter, which is assumed to be a perfect fluid and T^ v is 
the stress-energy tensor of the curvature fluid which represents an additional source 
term of purely geometrical origin [11] . By assuming /(i?) = exp(—R/A) we obtain 
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the field equations: 

G, v = -ke R /^l (5) 

(6) 



e -B/A 



\g^{A + R) + ^(R a R p + AR^) (g a „9 0v - g a pg^) 



In the case of the Friedmann-Lemaitre-Robertson- Walker (FLRW) metric, the above 
equations reduce to: 



(7) 



with, 



2 k H ■ R /A Ap \ 

d R H ■ 1 1 •• Ap . A 

2 a + 3-A R+ ^ R -A R -3^m {1 + Sw)+ 3 = °' 



/,'= .-<i[^+ij2 + A) . 

1 a a z 



where = - is the usual Hubble parameter and fc is the spatial curvature. The 
Bianchi identities applied to the total stress-energy tensor T™ T lead to the energy 
conservation equation for standard matter [2]: 

p + 3Hp(l + w) =0. (9) 
3. The vacuum case 

In the vacuum case (p — 0) equation |7]) and JHJ) can be written as a closed system of 
first order differential equations using the dimensionless variables, 

R R A k 

x = - — . v = z — K = — — - . (10) 

AH' y &H 2 6iJ 2 ' a 2 H 2 V ' 

It is clear that the variables y and z are a measure of the expansion normalized 

Ricci curvature and the cosmological constant respectively, K is the spatial curvature 

parameter of the Friedmann model, while a: is a measure of the time rate of change of 

Ricci curvature. The evolution equations for the variables (|10p are given by 

x' = 2z + 2K - 2 + x{l + x + y + K) , 
j/ =xz + 2y(2 + y + K) , 

z' = 2z(2 + y + K) , (11) 
K' = 2K(y + 1 + K) , 

where the prime represents the derivatives with respect to the time variable J\f = In a. 
This system is completed with the Friedmann constraint, 

l + K + y + z-x = Q, (12) 

which defines a hyperplane in the total phase space of the system. Consequently, all 
solutions of the dynamical system will be located in a non-compact submanifold of the 
phase space associated with (fTTj) . The time derivative of 172]) is nothing other than 
the Raychaudhuri equation. 
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3.1. Finite analysis 

The dimensionality of the state space of the system (fTTj) can be reduced by eliminating 
any one of the four variables using the constraint equation (Til?)) . If we choose to 
eliminate the x the dynamical equations become 

y' = y(4 + 2K + 2y + z) + z(l + K + z) , 

z' =2z(2 + K + y), (13) 
K' = 2K{l + K + y) , 

and the number of invariant submanifolds is maximized, making the analysis much 
easier. It is clear from the above equations that z = (corresponding to the (K — y) 
plane) is an invariant submanifold. Specifically, If we choose z ^ as an initial 
condition for our cosmological evolution, any general orbit can only approach z = 
asymptotically. This implies that no orbit crosses the (K — y) plane and consequently 
no global attractor can exist. 

Setting K' = 0, y' = 0, z' = 0, we obtain four fixed points. We can obtain exact 
cosmological solutions at these points using the equation, 

H = -(y + K + 2)H 2 . (14) 

In fact, at any fixed point, equation (| 14|) reduces to 

H = --H 2 , (15) 
a 

where 

a = (y» +K* + 2)- 1 , a^O (16) 

and the quantities X„ are meant to be calculated at the fixed point. Equation (TTJ} 
applies to both the matter and vacuum cases and describes a general power law 
evolution of the scale factor. In addition integrating with respect to time we obtain 

a = a (t-t ) a . (17) 

This means that by finding the value of a at a given fixed point, we can obtain the 
solutions associated with it using equation (TT4")) . 

In this way, points A v and B v (see Table [IJ are found to represent Milne 
and powerdaw evolutions respectively. However, by direct substitution into the 
cosmological equations it can be shown that these fixed points cannot be considered as 
physical points, which means that although we can choose initial conditions as close as 
we want to these points, the cosmology will never evolve according to these solutions. 

For the points C v and V v we have + + 2 = 0. In this case 

H = 0, (18) 

which implies 

a = a e 7(t -*° ) . (19) 

The value of the constant 7 can be obtained by direct substitution into equations ([7]) 
and (El) . For both C v and T> v we obtain 

(20) 



so they represent an exponential evolution. The contracting or expanding nature 
of this solution depends on the direction of approach of the orbits with respect to 
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the hypersurface y + K + 2 = 0. This hypersurface divides the phase space in two 
hypervolumes characterized by a contracting or expanding evolution. In particular, 
for y < —K — 2 the orbits describe a contracting universe, while for y > —K — 2 they 
represent an expanding one. 

The stability of the hyperbolic fixed points A v , B v and V v is obtained by using 
Hartman-Grobman theorem, while for the point C v , which is non- hyperbolic, we use 
the local center manifold theorem to find it's stability. A brief review of this approach 
in presented in the Appendix. 

In our case, using the transformation 

y = u — 2v — m , 

z = Am , (21) 
k = v 

the system (I13|) can be written in the diagonal form 

ii = - 4u + 12m 2 + 2mu + 2u 2 - imv - 2uv , (22) 
v = - 2v - 2mv + 2uv - 2v 2 , (23) 
m = — 2m 2 + 2mu — 2mv (24) 

near the non- hyperbolic fixed point C v . By substituting the expansions 

hi (to) = am 2 + bm^ + 0{m A ) , (25) 
h2(m) = cm 2 + dm 3 + 0(m 4 ) (26) 

into equations (j A.4|) and (|A.6[) and then solving for a, b, c and d, we obtain 

hl(m) = 3m 2 + ^m 3 + 0(m 4 ), h2(m) = 0(m 4 ) . (27) 

Substituting this result into equation (1A.4|) then yields 

m = -2m 2 + 0{m 3 ) (28) 

on the center manifold W c (0), near the point C v . This implies that the point C v is a 
saddle-node i.e. it behaves like a saddle or an attractor depending on the direction 
from which the orbit approaches. The local phase portrait in the neighborhood of C v 
is shown in Figure 1 . 

If one considers now the transformation (|21[) one realizes that m oc z, so that C v 
is an attractor for z > and a saddle for z < 0. This also clear from Figure [2] in which 
the invariant submanifold z — y is depicted. 

Finally, it is useful to derive an expression for the deceleration parameter q in 
terms of the dynamical variables: 

q = ~-l = -{y + K + l). (29) 

It follows that q > (y + K + 1) < 0. This condition is satisfied only for the point 
C v as expected by looking at the solution associated with this fixed point. In Figure 
[3] we give the location of the q = plane relative to the fixed points A v , C v and B v . 



3.2. Asymptotic analysis 

In this section we will determine the fixed points at infinity and study their stability. 
To simplify the asymptotic analysis we need to compactify the phase space using the 
Poincare method. Transforming to polar coordinates (r(N),6(N),<j>N)): 

z — > r cos 9, K — > r sin 9 cos <fr, y — > r sin 9 sin <f> 
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Table 1. Coordinates of the fixed points, eigenvalues, stability and solutions for 
exp(— ^)-gravity in vacuum. 



Point Coordinates(y,z,K) Eigenvalues Stability Solution 

A v [0,0,0] [2,4,4] repeller a = a a (t - t a ) 

B v [0,0,-1] [-2,2,2] Saddle a = a (t - t a )i 

C v [-2,0,0] [-4,-2,0] Saddle-node a = a e 7 (*-*°) 

V v [-2,1,0] [-M-^M] saddle a = a e^-^ 




Figure 1. The phase portrait for the system (12) in the neighbourhood of the 
fixed point C for exp (— i?/A)-gravity in vacuum. 



and substituting r — > yzr^, the regime r — > oo corresponds to R — > 1. Using this 
coordinate transformation and taking the limit i? — > oo, the system (13) can be 
written as 

-> j( 8cos(»sin(60 3 - sin(0) (-7sin(0) + sin(30)) 

+ 8cos(0) 2 sm(6»)(x) + 4 cos(0) sin(0) 2 sin^) ( x ) ) , (30) 

R0 , ^ -cos(0) 2 sin(0) (cos(fl)+sin(g) ( X )) ^ 
J? — 1 
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Figure 2. The phase space of the invariant submanifold z — y for exp (-R/ A)- 
gravity in vacuum. 



cos(cf>) cot(g) (x + cot(0)) 

^ S - ! ' ^ 

where x = cos((/>) + sin((/>) . Since equation (30) does not depend on the coordinate R, 
we can find the fixed points of the above system using equations (31) and (32) only. 
The results are shown in Table [2] 

Let us now derive the solution for the first fixed point, which corresponds to 
K — > — oo, z — * oo. In these limits the first equation of the system (13) reduces to 

K = 2K (33) 
and the equation for H becomes 

H = KH 2 . (34) 
Integrating equation (33) we obtain 

K = , . (35) 

Substituting K back into equation (34) and solving for (N — Noo), we obtain 

{N-N 00 ) = [c 1 ±~c (t-t )]§ . (36) 

The same procedure can be used to obtain solutions for the other asymptotic points. 
On the fixed circle 2^ we have three possible solutions corresponding to whether the 
sign of K and y are the same or different. 
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Figure 3. The location of the q = plane relative to the fixed points A v , C v and 
B v for exp (— i?,/A)-gravity in vacuum. 



The stability of the fixed points in Table[2]are found by expanding equations (31) 
and (32) up to second order and then applying the center manifold theorem to the 
resulting system (see Table [3]). 

For the fixed circle Too , the stability depends on the value of the angle 4>: 

-7r/4 < (j)< 7r/2 Stable, 

7r/2 < <p < 37r/4 Unstable, 

3?r/4 < (j) < 37r/2 Stable, 

37r/2 < <f> < 77r/4 Unstable. 

In the next section we will see how the introduction of matter modifies the picture we 
obtained in the vacuum case. 

4. The matter case 

In this case we can use the same dynamical variables we used for the vacuum case 
together with one additional variable D, that is related to the matter energy density: 

-A = JL A y- k n = Ap 

X AH' V 6# 2 ' Z 6# 2 ' a 2 H 2 ' 3ffV*/ A ' 1 ' 

The definition of the variables reveals that not all of the phase space corresponds to 

physical situations. This becomes clear if we divide D by z. We obtain 

®=2pexp(-j) , (38) 
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Table 2. Coordinates, eigenvalues and the Stability of the fixed points in the 
asymptotic regime for exp(— ^) gravity in vacuum. 



Point (9, (j)) Eigenvalues Coordinates Solution 

Aoo [^,0] [2,0] \K\,\z\^ooy^Q (N - N^) = [ Cl ± ^(t - t a )]i 

B v oo [f,-7r] [2,0] I^U^I^oot/^O (N-N 00 ) = [c 1 ±^(t-t )]i 

£ v oo [^L,0] — oo A"-* (JV - JV*,) = [ Cl ± f (t - t D )] 2 

[f^] [^:0] lyUzHooA-^O (AT - A^) = [ Cl ± - t )]f 

& 00 ^f'f] [^»°] li/I.N-^ooA-^O (JV-JV 00 ) = [ci±^(t-t )]' 

W„oo [f,|] [X, ] |y|,|z|->oo*-0 (7V-7V oc ) = [c 1 ±f (t-i )] 2 

Z„oo [■§,</>] [0, — x cos (0)] y ±oo X ^ =Foo a = a n exp[y(t — to)] 

y^+ooK^+oo (JV-JVoo) - [ Cl ±f (i-t )] 2 

y^-ootf^-oo (JV-JVoo) = [ci±^(t-t )]* 



Table 3. Coordinates of the asymptotic fixed points of exp (— _R/A)-gravity in 
vacuum and their stability. 



Point 


Coordinates(y,z,K) 


Stability 




[-tt/4,0] 


Unstable 




[tt/4, -tt] 


Stable 


f 


[-37r/4,-7r/2] 


Stable 


J oo 


[tt/4, -tt/2] 


Stable 


<?oc 


[-7r/4,7r/2] 


Stable 


Woo 


[37T/4, tt/2] 


Stable 



which has the same sign of p. This means that the sectors in the phase space for which 
the sign of D is different from the sign of z contain orbits in which standard matter 
violates the weak energy condition p > 0, and have to be discarded as not physical. 
As we will see this affect the sets of possible orbits for this model. 

Following the same procedure we used in the vacuum case, we obtain an 
autonomous system equivalent to the cosmological equations with non-zero matter 
density: 

x = 2 + 2z + 2K + x(l - x + y + K) - D(l + 3w) , 
y' = xz + 2y(2 + y + K) , 

z' =2z(2 + y + K), (39) 
K' = 2K(y + 1 + K), 
D' = D(l - 3w + 2y + 2K - x) , 
together with the constraint equation 

l + K + x + y-z-D = , (40) 
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where the prime again denotes the derivative with respect to the logarithmic time 
variable AT. 

4-1- Finite analysis 

The system (|39j) can be further simplified, by eliminating x using the constraint 
equation (j40|) : 

y' = y(4 + 2K + 2y + z) + z(l + K + D + +z) , 

K ' = 2K (1 + K + y) , (41) 

z' = 2z(2 + y + K) , 

D' = D(2 - 3w + 3K + D + 3y + z) . 

We have three invariant submanifolds: K = 0, z = and D — 0, so in this case also 
no global attractor can exist. Setting K 1 = 0, y' = 0, z' = and D 1 = we obtain 
seven fixed points. 

As in the vacuum case, we can use the coordinates of these fixed points and 
equation (|14p to find the behaviour of the scale factor at these points. In addition, 
the behaviour of the energy density p can be obtained from equation ([9|), which at a 
fixed point reads 

£ = -3(1 + (42) 
P t 

where a is defined by (|15[) . However, direct substitution in the cosmological equations 
reveals that all the fixed points correspond to vacuum states. 

Points A m and B m are found to represent Milne solutions while C m and V m 
represent a power law evolution. For points £ m , T m and Q m we find that H — 0, which 
means that these points represent Einstein-de Sitter solutions. The exact solutions at 
these fixed points are summarized in Table 2] 

As in the vacuum case, we use the Hartman-Grobman theorem together with the 
center manifold theorem to analyze the stability of all the fixed points. The results 
are shown in Table [5l 

Equation (|2"5|) . which relates the deceleration parameter to the dynamical 
variables generalizes in the matter case to the hyperplane, 

q = —^-l = -(y + K + l). (43) 
4-2. Asymptotic analysis 

We complete the analysis for the matter case by investigating the asymptotic behavior 
of the system (|4T|). In order to achieve this we compactify the phase space by 
transforming to 4-D polar coordinates. The transformation equations are 

D — > r cos 5, z^r sin 6 cos 9, x — > r sin S sin# cos</>, 

y — > r sinO sin 8 sin cj> , 

where r € [0, oof, S E [0, ir], 8 E [0, tt], and <fi E [0, 2tt]. We then transform the radial 
coordinate r — ► j^tr and in the limit R—*l, the system (40) reduces to 

R' -> cos((5) 3 + cos(£) cos(0) sin(<5) 2 sin(0) sin(0) + cos(J) 2 sin(<5) (cos(0) 

+ 3sin(0)x) +sin(5) 3 sin(6i)(cos(0)(2 + ^) (44) 
+ sin(0)(3 cos(6i) 2 + 2 sin(6>) 2 + ip) 
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Table 4. Coordinates of the fixed points, the eigenvalues, and solutions for 
exp (— _R/A)-gravity in the matter case. 



Point 


Coordinates(y,z,K,D) 


Eigenvalues 


Solution 




•A-m 


[0,0,0,0] 


[2-3io,2,4,4] 


a = a (t — 


to)? 




[0,0,-1,0] 


[2,2,-2,-(l + 3«;)] 


a = a (t — 


to) 


c 


[0,0,0,3u;-2] 


[3to-2,2,4,4] 


a = a a (t — 


to)? 


v m 


[0,0,-l,3w + l] 


[2,2,-2,(1 + 3^)] 


a = a a (t — 


to) 


f 


[-2,1,0,0] 


[ ^+ 3 ,^- 3 , 2, 3 2>w\ 


a = a e~ l(t 


-to) 


J~ rn 


[-2,0,0,0] 


[-2,-4,-310-4,0] 


a = a e~ l(t 


-to) 


Qm 


[-2,0,0,3™ + 4] 


[3^ + 4,-2,-4,0] 


a = a a e l(t 


-to) 



Table 5. Stability of the fixed points for exp (— i?A)-gravity in the matter case. 



Point 


w = 


< w < I 


w = \ 




Repeller 


Repeller 


Repeller 


B m 


Saddle 


Saddle 


Saddle 


Cm 


Saddle 


Saddle 


Saddle 


v m 


Saddle 


Saddle 


Saddle 




Saddle 


Saddle 


Saddle 


J~ m 


Saddle-node 


Saddle-node 


Saddle-node 


Qm 


Saddle-node 


Saddle-node 


Saddle-node 


Point 


3 <w< ! 


w = l 


| < W < 1 


Am 


Repeller 


Saddle-node 


Saddle 


B m 


Saddlee 


Saddle 


Saddle 


r 

*-*m 


Saddle 


Saddle-node 


Repeller 


v m 


Saddle 


Saddle 


Saddle 


f 


Saddle 


Saddle 


Saddle 


J m 


Saddle-node 


Saddle-node 


Saddle-node 


Qm 


Saddle-node 


Saddle-node 


Saddle-node 



RS' 



sin(<5) cos((5) 



8(R-1) 



8cos(S)(ip- 1) -sin((5)(cos(3l9) + 8cos(</>) sin(0) 

+ 8 sin(6>) 3 sm(<j>) + cos(0) (7 + 4 sin(6>) 2 (cos(20) 
- sin(2$))) 



R6 1 



cos(g) 2 
2(R-1) 



2 cos((5) sm(<f>) + sin((5)(2 cos(6>) sin(^) 
+ sin(0)(l - cos(2</>) + sin(2^)))l , 



(45) 



(46) 
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Table 6. Coordinates, eigenvalues, and the solutions for fixed points in the 
asymptotic regime for the exp(— R/A) gravity in matter case. 



Point (5, 9, tfi) Eigenvalues Solution 



A rv^ 


I 2 5 2 ' 2-1 


[1,0, OJ 




{AT AT \ _1_ 3c /j. 

(TV - iVoo) = [ci ± -^(i 
(N-nJ) = [a±%(t- 

\ / L J- 2 ^ 


— *o)J 3 

w / J 




[ 7T 3-7T 7T 1 

L 2 ' 2 ' 2-1 


[-1,0,0] 




(J m oo 


[" 7T 7T 7T] 

L 2 ' 4 ' 2-1 


[-£,0,0] 




(JV-JV 00 ) = [ci±^(t- 


T> m OQ 


I" 3-7T 7T 7T 1 

L 4 > 2 ' 2 J 


h£,o,o] 




(iV-AT 00 ) = [c 1 ±f(t- 


-a 2 


S m oo 


[ 2 ' 4 s 2 J 


[£,o,o] 




(JV-JV 00 ) = [c 1 ±%i(t 




J~ mOO 


L 2 ' 4 ' 2 -1 






( N — N ) — \r-i 4- f+ 


— t VI § 


M m O0 


[§,M 


[0,0, x ] 

y —* ±oo K — s 


■ ^Foo 


a = a ea;p[7(i - t )] 








j/ — > +oo if — i 


• +oo 


(N-N 00 ) = [c 1 ±f(t- 


-a 2 






y — ► — oo K — ! 


• — oo 


(N-N oc ) = [c 1 ±^(t 




O m oo 


[arccot(x),-^ A] 


[0,0,/(^)>0 

J/ — > ±00 if — i 


V0] 
■ ^Foo 


a = a exp[j{t - t )] 








J/ — > +00 if — i 


• +oo 


(N-N 00 ) = [c l ±%(t- 


-a 2 






y — ► — oo K — i 


• — oo 


{N-N 00 ) = [c 1 ±^(t 


"to)] 3 




;(</>) (cos((5) cot(0) - 


h cos(0) sin(<5)(sin 


w + 


cos(0) + cot(0))) 


(47) 






i?- 1 




! 



where (/J = cos(0) sin(6>) sin(</>). Notice that the first equation of the previous system 
does not depend on R, which means that the fixed points of this system can be 
determined by the angular equations alone. The solutions at the fixed points can then 
be obtained by following the same procedure we used in the vacuum case. All the 
results are summarized in Table [6] The stability of the first six fixed points are shown 
in TableH 

As before, we use the center manifold theorem to determine the stability of these 
fixed points. In this case, because the eigenvalues have a double zero, the coordinate 
which correspond to the non-zero eigenvalue are approximated by the function 

h(xi,x 2 ) =ax 2 + fex 2 x 2 + cx 2 + 0(|x| 3 ), (48) 

where Xi and x 2 are the coordinates that correspond to the zero eigenvalues (see the 
Appendix for details). The stability of the point A^ m oo depends on the angle cj): for 
— j > <fi > 2f , it is unstable, otherwise it is stable. Finally the point O m oo is unstable 
for all values of <fi. 

5. Discussion and Conclusions 

In this paper we have applied the dynamical system approach to the exponential 
gravity cosmological model, and found exact solutions together with their stability for 
both the vacuum and matter cases. 

In the vacuum case we found four finite critical points A v , B Vl C v and T> v , of 
which only two C v and T> v are found to be physical. These last points were found to 
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Table 7. Stability of the fixed points in non- vacuum exp (— R/A)-gravity. 



Point 


Stability 


A 


Stable 


-Boo 


Unstable 


c 


Unstable 




Unstable 


f 


Stable 


f 

J OO 


Unstable 



represent a solution whose nature depends on the parameter 7(A); for A > we can 
have either exponential expansion (7 > 0) or exponential contraction (7 < 0) and for 
A < the solution oscillates. 

From the stability point of view, the point C v , which lies in the invariant 
submanifold z — 0, is of particular interest because, since it is non-hyperbolic, it 
represent an attractor for z > and saddle for z < 0, while the other physical point 
V v is found to be a saddle. 

On the other hand, the solution connected with the non-physical points A v and 
B v are found to correspond to power law evolution and are also interesting because 
the orbits can approach arbitrarily close to them. 

In the asymptotic regime all the critical points represent solutions which have 
a maximum value for the scale factor, hence all models that evolve to one of the 
asymptotic future attractors will re-collapse. 

The invariant submainfold z = divides the phase space into two regions, z > 
and z < which correspond to A > and A < respectively. The fact that no orbit 
can cross the plane z = is then consistent with the fact that A is a fixed parameter 
for this model. 

In the vacuum case, we found that the region z < does not contain any finite 
critical point. However, in the plane z = 0, we have the physical point C v which 
represents a de-Sitter saddle and the non-physical points A v and B v are a repeller 
and saddle point respectively. Thus, the only attractors in the region z < are 
asymptotic, which means that all the models that begin their evolution in this region 
will re-collapse. It is also possible that one can choose initial conditions in such a way 
that both de-Sitter and power law phases are present in the evolution of the model. 

From a physical point of view, the region z > appears to be more interesting 
because the point C v represents a de-Sitter attractor for z > 0. It follows that there 
are two different possible solutions towards which models can evolve. 

Since the point T> v , which represents an unstable de-Sitter phase, lies in the region 
z > 0, a set of initial conditions exist for which orbits describe an intermediate de- 
Sitter phase (see figure [5]) . Furthermore, for models that evolve near the non-physical 
point B v an intermediate power law phase is also present. 

By looking at Figure 3 it is clear that the de-Sitter phases C v and V v are separated 
from the past attractor A v by the plane q = 0, therefore any model that starts near the 
past attractor A v and evolves toward the future de-Sitter attractor C v will cross the 
plane q = 0, indicating a transition from an accelerating evolution to a decelerating 
one. 
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Table 8. The sectors and the behaviour in each one for the exp(—R/A) 
gravity in matter case. Here +/— corresponds to a positive/negative values of 
the coordinates. 



y D z 


K Behaviour 


y D z 


K Behaviour 


+ - 


a 




- P 


+ + - 


a 




+ P 


+ - + 


a 


- - + 


1 


+ - 


+ a 


- - + 


+ 7 


+ + + 


a 


- + + 


e 


+ - + 


+ a 


- + + 


+ e 


+ + - 


+ a 


- + - 


- 8 


+ + + 


+ a 


- + - 


+ 5 



The introduction of matter into this model increases the dimensionality of the 
phase space, making it more difficult to visualize. On the bases of the relative stability 
of the fixed points and the invariant submanifold structure it is possible to catalog the 
possible evolutions of this model in 5 classes (see Table [8|). 

The first class (a) is characterized by the fact that that the cosmic histories 
evolve towards re-collapse. This class also contains cosmic histories which include 
an intermediate almost power law transient phase. For the second class (/3), cosmic 
histories evolve towards re-collapse, but in addition to a power law transient phase 
there can be also an oscillating one. The third class (7) contains two types of orbits. 
Depending on the initial conditions, the models will either evolve towards re-collapse 
or towards a de-Sitter type solution. During this evolution, a transient de Sitter phase 
can be present. A fourth class (e) also contains two types of orbits. The universe 
can either re-collapse or evolve to a de-Sitter type model. In this case there are 
transient phases that include two different unstable de-Sitter evolutions . The final 
class (S) contains models that either re-collapse or end up at an oscillating solution. 
The possible transient phase in this class are multiple oscillations and/or almost power 
law behaviour. As it can be seen from Tabled the condition D/z > coming from the 
weak energy condition, required for standard matter, allows one to exclude completely 
the classes 7 and S. 

In conclusion, exp (—7^) gravity has a very rich structure that includes a series of 
diverse cosmological histories. Particularly important are the ones including multiple 
de Sitter phases because they could provide us with natural models describing the early 
and late time acceleration of the Universe. Unfortunately, as is clear from Figure 
this scenario does not include a decelerated expansion phase between these two de- 
Sitter phases, so unless some other mechanism is taken into account in these cosmic 
histories will not admit a standard structure formation scenario. 

Appendix A. 

In section 3.1 we used the center manifold theorem to analyze the stability of the 
non-hyperbolic fixed point C v in the vacuum case, here we will give a brief review of 
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this approach [12] . Consider the nonlinear system (bold letters represent vectors) 

u = /(u). (A.l) 

For simplicity we shall assume that the origin is a non-hyperbolic fixed point for this 
system (this assumption does not affect the generality of our treatment because it is 
always possible to change the coordinates to make the fixed point the origin of the 
new coordinate system). If / € C 1 (i5) and /(0) = 0, then this system can be written 
in the diagonal form 

u=Ju + T(u), (A.2) 

where J = Df(0) — diag[Z, P, N], the square matrices Z,P,N have r eigenvalues of 
zero real part, p eigenvalues of positive real part and n eigenvalues of negative real 
part respectively and T(u) = /(u) - Ju, where T E C^E), T(0) = and DT(0) = 0. 
The system (48) can be divided into three subsystems 

x = Zx + F(x,y,z), 

y = Py + G(x,y,z), (A.3) 
z = Az + ff(x,y,z), 

where (x,y,z) e R r x BP x R n , F(0) = G(0) = H(0) = 0, and L>F(0) = DG(0) = 
DH(0) = 0. If (F, G, H) G C m {E) with m > 1, it follows from the local center 
manifold theory that there exist a z-dimensional invariant center manifold Wf ocal (0) 
defined by 

W^ ocq; (0) = {(x,y,z) e R r xPx R n \y = hl(x),z = h2(x) for|z| < 5} , 

for some S > 0, where (hl,h2) € C r (N s (0)), hl(0) = h2(0) = 0, Dhl(0) = 
Dh2(0) = 0, and they satisfy 

Dhl(x)[Zx + F(x, hl(x), h2(x))] - Phl(x) - G(x, hl(x), h2(x)) = , 
L>h2(x)[Zx + F(x, hl(x), h2(x))] - A^h2(x) - H(x, hl(x), h2(x)) = . 

In the neighborhood of a non-hyperbolic fixed point the qualitative behavior of the 
system (48) is equivalent to the qualitative behavior of the reduced system 

x = Zx + F(x,hl(x),h2(x)). (A.4) 

The functions hl(x), h2(x) can be approximated by substituting the series expansion 
of their components into equations (51) and (52). In the case when we have a double 
zero eigenvalues the coordinate which correspond to the non-zero eigenvalue can be 
approximated by the function, 

h(xi,x 2 ) -ax2 + 6 x 2 x 2 + cx2 + 0(|x| 3 ), (A.5) 

where Xi and x 2 are the coordinates which correspond to the zero eigenvalues. In 
general the flow on the center manifold near the fixed point takes the form 

x = ax r + ..., (A. 6) 

If r > 2 and a r ^ 0, then for r even we have saddle-node at the fixed point, for r odd 
and a r > we have unstable node and for r odd and a r < we have a topological 
saddle. 
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